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ABSTRACT 

We present results from two observations of the wind-accreting X-ray pulsar 4U 1907+09 using the Suzaku 
observatory. The broadband time-averaged spectrum allows us to examine the continuum emission of the 
source and the cyclotron resonance scattering feature at ~19keV. Additionally, using the narrow CCD response 
of Suzaku near 6 keV allows us to study in detail the Fe K bandpass and to quantify the Fe K/3 line for this 
source for the first time. The source is absorbed by fully-covering material along the line of sight with a column 
density of Nh ~2 x 10 22 cm" 2 , consistent with a wind accreting geometry, and a high Fe abundance (~3-4x 
solar). Time and phase-resolved analyses allow us to study variations in the source spectrum. In particular, dips 
found in the 2006 observation which are consistent with earlier observations occur in the hard X-ray bandpass, 
implying a variation of the whole continuum rather than occultation by intervening material, while a dip near 
the end of the 2007 observation occurs mainly in the lower energies implying an increase in Nh along the line 
of sight, perhaps indicating dumpiness in the stellar wind. 

Subject headings: X-rays: stars — X-rays: binaries — stars: pulsars: individual (4U 1907+09) — stars: 
magnetic fields 



1. INTRODUCTION 

Discovered by Giacconi et al. (1971), the X-ray pulsar 
4U 1907+09 is a wind-accreting high mass X-ray binary 
(HMXB) with a highly reddened companion star of magni- 
tude mboi = 16.37 mag, luminosity Lboi = 5 x 10 5 L Q and a 
mass loss rate of M = 7 x 10" 6 M yr" 1 (Cox et al. 2005). The 
neutron star is in an eccentric orbit with a period of 8.3753 d. 
Orbital parameters can be found in in't Zand et al. (1998). 

From 1983, when Makishima et al. (1984) found the source 
to have pulsations, through 1998 the pulse period had a steady 
spin down rate of f pu i se = +0.225syr _1 , going from 437.5s 
to 440.3 s (in 't Zand et al. 1998). Baykal et al. (2006) re- 
ported a deviation from the spin down rate, measuring a much 
lower value between 1998 and 2003, which in 2002 was half 
the long-term value (~ 0. 1 15 syr" 1 ). Using INTEGRAL data, 
Fritz et al. (2007) demonstrated that between 2004 and 2005 
a complete torque reversal occurred with the maximum pe- 
riod of ~441 .3 s being reached in April 2004 followed by spin 
up behavior with P pu i se = -0.158 ± 0.007 syr" 1 . Inam et al. 
(2009), using RXTE -PC A data between 2007 and 2008, re- 
ported a second torque reversal and a new spin down rate of 
+0.220syr _1 , consistent with the rate before 1998. 

Using data from Ginga observations of 4U 1907+09 Mak- 
ishima et al. (1992, 1999) reported a cyclotron resonance 
scattering feature (CRSF) at ~19keV. CRSFs or "cyclotron 
lines" are generated by photons scattering on electrons in the 
strong magnetic field near the surface of a neutron star. Pho- 
tons from the accretion column above the magnetic poles that 
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are near the energies that separate the Landau levels of the 
electrons are scattered out of the line of sight, the effect re- 
sembling an absorption line (see, e.g., Nagel 1981; Schonherr 
et al. 2007 for details). These features can thus be used to di- 
rectly measure the magnetic field strength of the pulsar. Using 
BeppoSAX observations Cusumano et al. (1998) later con- 
firmed the CRSF for 4U 1907+09 at -19keV and calculated 
a strong surface magnetic field of 2.1 x 10 12 G. Cusumano et 
al. also reported a harmonic at ~39 keV and an iron emission 
line near 6.4 keV with an equivalent width of about 60 eV. 

4U 1907+09 was discovered to have an orbital period of 
8.38 days by Marshall & Rickets (1980) and to display peri- 
odic flaring behavior over the orbit as well as occasional dip- 
ping behavior. The flaring was originally believed to be due 
to a Be type companion star, where the flares corresponded to 
the neutron star passing through the companion's outer atmo- 
sphere. Van Kerkwijk et al. (1989), however, optically iden- 
tified the companion star 4U 1907+097 as a supergiant, con- 
tradicting the earlier assumption. This was later confirmed 
by Cox et al. (2005) who identified the companion as a type 
08/09 supergiant with a dense stellar wind, leading to the 
current belief that the pulsar is entirely wind-accreting. 

The dipping behavior has been observed, notably by the 
Rossi X-Ray Timing Explorer (RXTE), as reductions in in- 
tensity for periods on the scale of an hour, occasionally drop- 
ping completely below detection levels (in't Zand et al. 1997). 
Due to the absence of strong variations in the column den- 
sity, the lack of correlation with the orbital phase, and various 
spectral features, in't Zand et al. have suggested that these 
dips are due to the cessation of mass accretion rather than oc- 
cultation of intervening materal. 

In this paper we present data from the Japanese X-ray ob- 
servatory Suzaku which observed the source twice, in 2006 
and 2007. In §2 we describe the instruments and the data re- 
duction process. In §3 we describe our methods and results for 
time-averaged spectral fitting, including details of the 6.4 keV 
iron emission line and the CRSF at 19keV. §4 shows the re- 
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suits of pulse-period determination and phase-resolved spec- 
tral fitting. In §5 we detail time-resolved spectral fitting in 
order to explore the dipping behavior seen in both observa- 
tions and the flaring seen in the 2007 observation. Sections 
§6 and §7 are devoted to a detailed discussion and a summary 
of our conclusions respectively. 

2. DATA REDUCTION AND ANALYSIS 

Suzaku observed 4U 1907+09 with the X-ray Imaging 
Spectrometer (XIS; Mitsuda et al. 2007) and the Hard X-ray 
Detector (HXD; Takahashi et al. 2007) beginning 2006 May 2 
at 06:11 UT (MJD 53857.3) for 123 ks and again 2007 April 
19 starting at 10:03 UT (MJD 54209.4) for 158 ks. The obser- 
vation ID's are 401057010 and 402067010 respectively. The 
HXD gathered data with both its detectors, the PIN diodes and 
the GSO scintillators, however we did not use the GSO data 
because of the faintness of the source in the GSO band rela- 
tive to the non-X-ray background. The XIS has 4 CCDs: 0-3 
with 0, 2 and 3 front-illuminated and XIS 1 back-illuminated. 
In November of 2006 XIS2 was damaged and taken off-line, 
therefore data in the 2007 observation were taken with only 
the remaining three XISs. Two corners of each XIS CCD have 
an 55 Fe calibration source which can be used to calibrate the 
gain and test the spectral resolution of data taken using this 
instrument (See the Suzaku Data Reduction Guide 5 for de- 
tails). 

During the 2006 observation the XIS was set to 
1/4-window mode to avoid saturating the CCD should the 
flux during the source's flares have proven to be too high. 
However, this did not prove to be the case and for the 2007 ob- 
servation the mode was switched to full-window. The point- 
ing was also changed from XIS nominal pointing in 2006 to 
HXD nominal, both to reduce the flux in the XIS by a small 
fraction (10%) and to increase the sensitivity above 12keV 
for analysis of the CRSF. 

Data were processed with version 2.0.6.13 of the Suzaku 
pipeline and screened to exclude data within 436 s of passing 
through the South Atlantic Anomaly and within Earth eleva- 
tion angles of less than 5° or Bright Earth angles of less than 
20°. All extractions were done using HEASOFT v.6.6.2. 

2.1. XIS Reduction 

After screening, the good exposure time per XIS was 
58.4 ks for 2006 and 80.6 ks for 2007. The XIS CCDs were 
in the 3x3 and 5x5 editing modes which were cleaned 
and added to create image files for each XIS. Data from 
each CCD were analyzed separately. From these we extraced 
source and background lightcurves and spectra and used XIS- 
RMFGEN and XISSIMARFGEN to create the response ma- 
trix (RMF) and ancillary response (ARE) files. There was 
no charge injection in either observation. Channels were 
grouped by 5 channels below 1.5 keV and 10 channels from 
1.5-10keV. Data below 1 keV and above lOkeV were ig- 
nored for all XIS in both observations, due to complete ab- 
sorption below 1 keV and instrument noise above lOkeV. In 
addition, data between 1.8 and 2.7 keV and between 1.65 and 
2.3 keV were ignored for the 2006 and 2007 observations re- 
spectively, due to poorly understood calibration for the Si K 
complex and Au M edge arising from the detector and mirror 
system. 

Fitting the 55 Fe calibration sources for 2007 in XSPEC 
v.ll.3.2ag (Arnaud et al. 1996) with a model comprised 

5 http://heasarc.gsfc.nasa.gov/docs/suzaku/analysis/abc/abc.html 
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FIG. 1. — Time-averaged spectral fitting for 2006. Panel a) shows the data, 
all four XISs and the PIN. Panel b) the data/model residuals for an absorbed 
continuum model with FDCut. Panel c) shows the data/model residuals for 
our best-fit model including the iron lines, lOkeV feature and the CRSF. 

of three Gaussian components (Mn Ka\, Kct2 and K/3) 
yielded the following results for the Mn Ka\ line energy 
(expected value of 5.899 keV): 5.829 ± 0.004 keV for XISO, 
5.886 ± 0.006 keV for XIS 1 and 5.920 ± 0.007 keV for XIS3, 
showing that the energy calibration is good to 20 eV in XIS 1 
and XIS3 and to 70 eV in XISO. No calibration spectra were 
available for the 2006 data which was in 1 /4-window mode 
and therefore did not see the corners of the CCDs. 

2.2. PIN Reduction 

The Suzaku HXD/PIN is a non-imaging instrument with 
a 34' square (FWHM) field of view. The HXD instru- 
ment team provide non-X-ray background (NXB) event 
files using the calibrated GSO data for the particle back- 
ground monitor ("background D" or "tuned background" with 
METHOD=LCFITDT). This yields instrument backgrounds 
with <1.5% systematic uncertainty at the 1<t level (Fukuzawa 
et al. 2009). As suggested in the Suzaku ABC Guide, the 
Cosmic X-ray Background was simulated in XSPEC v.l 1 us- 
ing the form of Boldt (1987). 

Net spectra were extracted and deadtime-corrected for a net 
exposure time of 32.0 ks for 2006 and 64.9 ks for 2007. We 
excluded PIN data below 13 keV due to thermal noise. In both 
observations the source flux falls to levels well below a few 
percent of the background at ~40keV for our time-averaged 
spectra and therefore we only use data between 1 3 and 40 ke V. 
Channels were grouped with a minimum of 50 counts per bin 
to allow x 2 fitting. 

3. TIME- AVERAGED SPECTRAL ANALYSIS 

All spectral fitting in this paper used XSPEC v.ll.3.2ag 
(Arnaud et al. 1996). Unless otherwise stated, all uncertain- 
ties below correspond to Ay 2 = 2.71 (90% confidence level 
for one parameter). For the time-averaged fits the normaliza- 
tion for XISO was fixed while the cross-normalizations were 
left free for XIS 1 , 2 & 3 (2006) and for XIS 1 & 3 (2007); best- 
fit values of x 2 /dof were usually close to unity. The PIN/XIS0 
cross-normalization was set to 1.16 (2006, XlS-nominal) or 
1.18 (2007, HXD-nominal pointing; Maeda et al. 2008). Al- 
lowing the cross-normalization to be free resulted in unreal- 
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FIG. 2. — Time-averaged spectral fitting for 2007. Panel a) shows the data, 
three XISs and the PIN. Panel b) the data/model residuals for an absorbed 
continuum model with FDCut. Panel c) shows the data/model residuals for 
our best-fit model including the iron lines, lOkeV feature and the CRSF. 
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FIG. 3. — Data/model residuals focusing on the iron K bandpass for 2006 
when the continuum is modeled as in Fig. 1 b. 
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FIG. 4. — Data/model residuals focusing on the iron K bandpass for 2007 
when the continuum is modeled as in Fig. 2b. 



istic values due to a broad feature around 8-15 keV. We fit 
the XIS+PIN data in the range 1-40 keV. The 2006 data are 
shown in Fig. la; the 2007 data are shown in Fig. 2a. 

Our initial continuum model consisted of a power-law con- 
tinuum emission component with a high-energy rollover, plus 
full-covering absorption by a cold column of gas along the 
line of sight. For both the 2006 and 2007 observations, we 
considered three forms for the power-law continuum rollover: 
a power-law with a smooth Fermi-Dirac cutoff "FDCut" 
(Tanaka 1986), a power-law with an abrupt exponential cutoff, 
"HighECut" (White, Swank & Holt 1983), and the "NPEX" 
model consisting of negative and positive power-laws with a 
common exponential cutoff (Mihara 1995). Qualitatively, all 



models gave similar overall results for each observation. For 
the purposes of this paper we will be using the FDCut model 
to illustrate how spectral fitting proceeded and to fit the time- 
and phase-resolved spectra. We include the other two models 
as an alternate phenomenological descriptions of the contin- 
uum and for ease of comparison with other observations. 

For the absorber, we assumed the elemental abundances 
of Wilms, Allen & McCray (2000) and the photo-absorption 
cross sections of Verner et al. (1996); solar abundances were 
initially assumed. Data/model ratios to this best-fit initial FD- 
Cut model are shown in Figs.QJ) and|2]3 for 2006 and 2007, re- 
spectively. Initial values for x 2 /dof were 2689.41/996 (2006) 
and 3820.6 1/784 (2007). 

Figs.|3]and|4]show the data/model ratios in the Fe K band- 
pass. Both spectra showed narrow, positive data/model resid- 
uals at 6.4 keV, evidence for the Fe Ka emission line. We 
modeled this line using a Gaussian component. In the 2007 
fit we allowed the line energy to be free for all three XIS's 
based on the differences in gain observed in the calibration 
source. This allowed for a much better fit of the Ka line in 
the 2007 data. It was not necessary to do this in the 2006 data. 

Strong data/model ratios near ~15— 25 keV (Figs.QJ) &HJ)) 
indicated the presence of the CRSF, which we modeled using 
a Gaussian component in absorption. The actual shape of the 
CRSF is more complicated, but for the purposes of fitting, a 
Gaussian is adequate to model the feature (Schtinherr et al., 
2007). In our best-fit model the Gaussian centroid energy was 
18.6__^keV (2006) and 19.1 ±0.2 keV (2007), consistent with 
previous results by Cusumano et al. (1998) and Makishima et 
al. (1992, 1999); note that Makishima et al. used a Lorentzian 
line profile while we and Cusumano et al. used a Gaussian, 
however the two profiles have similar energy centroids de- 
spite their different shapes. We were not able to search for the 
second harmonic at ^39 keV since the PIN spectrum becomes 
background dominated in that region. 

We modeled positive ratios above 8 keV in the XIS and 
below ~15 keV in the PIN with a broad emission feature 
peaking around 10-13 keV. The "10 keV" bump feature has 
been observed previously in this source by Mihara (1995) 
and Coburn et al. (2002) as well as in several other cyclotron 
line sources, e.g., with RXTE in MXB 0656-072 (McBride 
2006), Her X-l, 4U 1626-67, and 4U 1538-52 (Coburn et 
al. 2002), and Cen X-3 (Suchy et al. 2007); and with Ginga 
in 4U 1538-52, and V 0331+53 (Mihara 1995). Though the 
physical nature of this feature is still not clear, it cannot stem 
from calibration issues since it has been seen with multiple 
satellites. It is believed to be part of the overall continuum 
and including this feature as a broad Gaussian emission com- 
ponent in the model is required for a good fit. x 2 /dof fell from 
1660.25/990 to 1125.39/987 (2006) and from 1743.56/776 
to 1315.54/773 (2007). The best-fit energy centroid for this 
component falls within or close to the gap in energy coverage 
between the XIS and PIN and the effective areas of the XIS 
and PIN are low in this region. The energy, width and inten- 
sity of this component as derived with Suzaku may therefore 
be affected by systematic uncertainties which are difficult to 
quantify. 

Next, we noticed narrow positive data/model ratios near 7. 1 
keV, seen in all four XlSes in 2006 and primarily in XIS 1 and 
3 in 2007. We identify this feature as Fe K/3 emission; detec- 
tion of this line in 4U 1907+09 is reported here for the first 
time. Modeling it with a Gaussian component, with the en- 
ergy left as a free parameter and the width of the K(3 line kept 
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FIG. 5. — Pulse period evolution from 1983 to 2008, using data from Baykal 
et al. (2001, 2006), Fritz et al. (2006), Inam et al. (2009) and this work (indi- 
cated by diamonds). 

tied to that of the Ka line, x 2 /dof fell to 1090.1 1/986 (2006) 
or 1290.08/770 (2007). In our best-fit model, it was signif- 
icant at the 99.994% (2006) or 99.998% (2007) confidence 
level to include this component according to an F-test (which 
formally cannot be used in this manner due to a boundary con- 
dition violation as per Protassov et al. (2002), but can give a 
rough approximation of the significance). The best-fit value 
of the intensity was 13 ±5% that of the Ka line, consistent 
with an origin in neutral or at most moderately-ionized (< 
XIII) species of Fe (Palmeri et al. 2003). For the remainder 
of this paper, the K/3 normalization is fixed to 13% that of the 
Ka line. 

Finally, we noticed absorption-like data/model ratios be- 
tween 7-8 keV. We fit this edge-like feature by allowing the 
Fe abundance of the absorber relative to solar, Zp e , to vary. 
As shown in Table Q] values of Zp e near 3^4 were typically 
obtained. We thus obtained our best-fit model, with x 2 /dof 
falling to 1040.19/985 (2006) or 1239.02/766 (2007). This 
edge is consistent with our values for the Fe Ka and K/3 lines, 
all of which have significant physical implications as will be 
discussed in §6.3. 

Figs, lc and 2c show the data/model ratios for our best- 
fit FDCut models, with all the above-mentioned components 
included. Best-fit model parameters are listed in Table Q]for 
the FDCut, HighECut, and NPEX fits to the 2006 and 2007 
spectra. Flux values for the two observations were 1 .5 and 2.8 
x 10~ 10 ergs cm" 2 s" 1 . The fits are robust to ±2% variations in 
the PIN background level (Fukazawa et al. 2009). 



4. PHASE RESOLVED SPECTROSCOPY 

4.1. Pulse Period Determination 

To measure the pulse period we performed epoch folding 
(Leahy et al. 1983) on the PIN data for each observation, 
excepting the gaps due to the satellite orbit and only using 
the first ~ 50 ks for the 2006 observation to avoid the pro- 
nounced dip beginning midway through that observation. The 
2006 observation pulse period was found to be 441 .09±0.05 s 
while the 2007 observation pulse period was found to be 
441.03±0.03 s. Uncertainties are due to strong pulse-to-pulse 
variations as well as the brevity of our observations compared 
to the length of the pulse period. Nevertheless, these num- 
bers are consistent with the trends in the pulse period seen by 
both Fritz et al. (2006) and inam et al. (2009). The long-term 




FIG. 6. — Energy resolved pulse profiles for the co-added XIS and PIN data 
from 2006. Flux contributions are listed for each energy band, representing 
the percentage of the total flux contributed by that energy range. 
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FIG. 7. — Energy resolved pulse profiles for the co-added XIS and PIN data 
from 2007. Flux contributions are listed for each energy band, representing 
the percentage of the total flux contributed by that energy range. 

evolution of the pulse period is shown in Fig. [5] with the two 
Suzaku observations falling in the midst of the second torque 
reversal, restoring the source to its previous long-term spin- 
down rate (inam et al. 2009). Physical implications of this 
event will be discussed in §6. 

Using these periods we generated pulse profiles for each 
observation. The profiles are shown in Figs. [6] and [7] exhibit- 
ing a double-peak with the initial peak being much weaker in 
the higher energy bands. 
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TABLE 1 

Time-Averaged Spectral Fit Parameters 



Model FDCut 



Parameter 

2-10 keV Flux (ergs cm -2 s _1 ) 
N H (10 22 cm- 2 ) 
Fe Abundance (rel. to Solar) 
Photon Index 
Power Law Norm 1 
Cutoff Energy (keV) 
Folding Energy (keV) 
NPEX a 1 
NPEX al 

NPEX Norm 1 (1(T 2 ) 
NPEX Norm 2(10- 3 ) 
NPEX Temperature (keV) 
CRSF Energy (keV) 
CRSF a (keV) 
CRSF Depth 

10 keV Bump Energy (keV) 
lOkeV Bumpcr(keV) 
10 keV Bump Intensity 2 
Fe Ka Line Energy (keV) 
Fe Ka Line a (eV) 
Fe Ka Line Intensity 2 
Fe Ka Line Equivalent Width (eV) 
Fe K/3 Line Energy (keV) 
X 2 /dof 
Reduced \ 2 

Units are (phkeV- 1 cm" 2 s _I at 1 keV) 
2 Units are (phcm~ 2 s" 1 ) 

4.2. Phase Resolved Spectral Analysis 

We divided the pulse profiles into 6 regions as shown in 
the first sections of Figs.[H](2006) and|9](2007) and extracted 
spectra for each region separately. Each bin had an average 
good exposure time per XIS of 9.5ks for 2006 and 13.2ks 
for 2007. The number of bins was chosen as a trade-off be- 
tween achieving good time or phase resolution and having ad- 
equate signal/noise within each sub-spectrum. We applied the 
best-fit time-averaged FDCut models, allowing the following 
parameters to be free: the power-law photon index (T), the 
power-law normalization, the cutoff energy, the CRSF energy 
and depth, the Fe line intensity and the 10 keV bump energy 
and depth. Other parameters were kept fixed in the fits at their 
time-averaged values; for some, such as the e-folding energy, 
there was no evidence of variability when these parameters 
were left free, thus we kept these parameters fixed to mini- 
mize any potential systematic effects on other free parame- 
ters. The iron line energies and widths were also kept fixed 
since we do not expect them to vary with phase. The results 
for selected model parameters are plotted in Figs. [8] (2006) 
and [9] (2007). We caution, however, that with a low number 
of pulse bins, evidence for variation in these parameters as 
well as correlation with flux is not highly robust. 

For both sets of data we find that the folding energy does 



HighECut NPEX 



not exhibit strong evidence of variation over the pulse phase. 
The cutoff energy and T do seem to vary over the phase, with 
both parameters reaching high values in bins and 5 and min- 
imum values near bins ^2-4. The Fe line depth shows some 
variation over the pulse with minimum values in the largest 
dip. There also seem to be small variations in the CRSF pa- 
rameters, all of which reach their largest values during the first 
pulse peak. However, given the low number of counts above 
lOkeV during the first peak and the influence of the parame- 
ters of the 10 keV feature produce a lack of robustness in this 
result. The fit results also suggest that the properties of the 10 
keV feature may vary as a function of pulse phase, but as we 
cautioned earlier, there may be systematic effects hampering 
our ability to measure the properties of this feature accurately. 

5. TIME RESOLVED SPECTROSCOPY 

For a preliminary exploration of spectral variability as a 
function of time, we present count rate light curves binned on 
the timescale of the satellite orbit in Figs. [10] (2006) andfTTI 
(2007). Softness ratios are also shown, including ( 1— 3)/(3— 
5)keV, sensitive to variations in Nn, and (3-5)/(5-10)keV, 
sensitive to variations in T. Additionally, Figs. [12] and [13] 
show hardness ratios plotted against each other for each inter- 
val. The observational data points in these color-color plots 
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FIG. 8. — Results for phase resolved spectral fitting for 2006. Pulse profile 
is for XIS0. For units see Table 2. A star symbol (*) indicates that the energy 
was unconstrained because the component was not needed for good fit. The 
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12.5 keV in order to avoid degeneracy between these two components. 
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was unconstrained because the component was not needed for good fit. 
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were obtained adding data from XIS0 and XIS3, while the 
grids were created by deriving colors from simulated spectra. 
The latter are based on the spectral parameters of the point 
with the highest value of 1.5-3keV over 3-10keV (exclud- 
ing points from the dips and flares) and allowing for an ad- 
ditional absorption component, partially covering the X-ray 
source, i.e., the "pcfabs" component in XSPEC. Using the re- 
sponse matrices for XIS0 and XIS3, hardness ratios for a grid 
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of values of the covering fraction and the additional Nu were 
calculated and then interpolated to create the solid and dotted 
lines in Figs. [T2l and [T3l 

For the 2006 observation, the softness ratios vary over time, 
beginning and ending hard, while softening in the middle, 
suggesting that Nh is varying over the observation. The pro- 
longed dip in 2006 between 60 ks and lOOks shows an in- 
teresting behavior, beginning very soft but then hardening 
quickly. These points are shown in Fig. [12] as squares and 
clearly do not lie along the partial covering grid, indicating 
that this dip is not caused by a dramatic increase in Nu- To- 
gether with the softness at the beginning of the dip, this im- 
plies a variation in T and that the cause of the dip is not ab- 
sorption but is primarily associated with the continuum com- 
ponent. For the 2007 observation, the (l-3)/(3-5)keV soft- 
ness ratio light curve exhibits a dip during the last ~40 ks of 
the observation, suggesting an increase in Nu, possibly caused 
by, for example, a relatively dense clump of material travers- 
ing the line of sight. These dip points (indicated by squares 
in Fig. 1 1 31> are consistent with a high density partial covering 
material in the line of sight. The (3-5)/(5-10)keV softness 
ratio light curve is roughly constant, suggesting that V does 
not strongly vary, except possibly towards the final ^20 ks. 

We divided the 2006 and 2007 data into 6 and 14 bins 
respectively, trying to ensure enough counts in each bin to 
make spectral fitting possible, resulting in each bin having an 
average good exposure time of roughly 1 1 ks (2006) or 5 ks 
(2007). The significant difference in binning for the two ob- 
servations stems from the extended dip in 2006 (see the top 
panel of Fig. [T4"l i which caused the average flux to be ^2.5 
times fainter than 2007. We applied the best-fit time-averaged 
FDCut models, keeping the iron line energies and widths, the 
lOkeV feature energy and width, and the e-folding energy 
fixed at the time-averaged value. The light curves and se- 
lected model component parameters are plotted in Figs. [14] 
(2006) andQ~5](2007). 

In both observations we notice that the Fe line depth is 
strongly correlated with the flux. Correlation coefficients 
were found to be 0.86 for 2006 and 0.95 for 2007 with null- 
hypothese probabilities of 97% and > 99.99% respectively. 
The implications for the geometry of the Fe line-producing 
material around the source will be discussed in §6. 

For the 2006 observation we do see some variation in 
Nh and T, particularly during the dip between 60 and lOOks. 
2D contour plots of Nh versus T calculated with the cutoff 
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energy left as a free parameter for the lowest flux subspec- 
trum show no strong correlation between the two parameters; 
additionally Monte Carlo simulations show a distribution in 
parameter values only on scales similar to the calculated error 
bars which are very small compared to the observed parameter 
variations. Standard deviations for A^h and T were 0.02 cm" 2 
and 0.015 respectively. In addition, we investigated possible 
correlations between T and the cutoff energy in the low flux 
states but again 2D contour plots and Monte Carlo simulations 
indicated that the degree of correlation was much smaller than 
the observed variation. 

For the 2007 observation there is a doubling in A^h during 
the final 2-3 spectra, from — 2 x 10 22 to ~ 4 x 10 22 cm" 2 . 
Through most of the observation T is roughly constant, but 
then it too changes in the final 3 spectra, flattening signifi- 
cantly where the flux reaches its lowest values. Contour plots 
of T as a function of A^h for the final three spectra yielded 
no evidence for any strong degeneracy between these two pa- 
rameters from one time bin to the next (Fig.fToTl. It is possible 
however, that the dip in T may actually be a result of not us- 
ing a partial covering model for the absorber, since we see in 
Fig. [13] that the dip is consistent with partial covering. The 
cutoff energy also decreased during this time, but we caution 
that there may be systematic effects associated with the lowest 
flux levels. There is no strong evidence for variation in the pa- 
rameters of the CRSF or the lOkeV bump energy and width, 
however the lOkeV bump normalization appears to track the 
continuum flux, consistent with this feature being a part of the 
emission continuum we do not yet fully understand. 

6. DISCUSSION 
6.1. Properties of the CRSF 

Observations of CRSFs allow direct measurement of the 
neutron star's magnetic fields, as the line energy E cyc oc B/( 1 + 
z), where z is the gravitational redshift of the line-forming re- 
gion. Their properties (depths, profiles, etc.) depend on the 
pulsar X-ray emission processes and the behavior of matter 
in strong magnetic fields (Meszaros 1985). CRSFs thus carry 
a wealth of information about the environment in which they 
are produced. For 4U 1907+09 we measure the energy cen- 
troid of the CRSF to be consistent with previous measure- 
ments, e.g., that of Cusumano et al. (1998) who calculated a 
magnetic field of ~ 2.1 x 10 I2 G using a value of (1 +z)~ l = 
0.76. The depth and width are also consistent with previous 
measurements. 

In recent years it has been discovered that in most accret- 
ing neutron stars the cyclotron line energy depends on the 
source flux (Mihara 2004, 2007; Nakajima 2006). Promi- 
nent recent examples for such a variability are V0332+53 
and 4U 0115+63, where the cyclotron line energy decreases 
with flux (Mowlavi 2006; Tsygankov 2006; Nakajima 2006), 
and Her X-l, where the line energy increases with flux 
(Staubert et al. 2007). The most prominent counter-example 
is A0535+262, where the line does not show any flux depen- 
dence (Caballero et al. 2007). The line variability can be un- 
derstood by considering that for many accreting neutron stars, 
even though the overall source luminosity is only at a few per- 
cent of Eddington, the radiation field produced in the accre- 
tion column is close to or above the local Eddington limit. In 
this case, the accretion stream is rapidly decelerated by radi- 
ation pressure from the magnetic poles leading to a radiative 
shock, close to which the cyclotron lines are formed (Becker 
& Wolff 2005, 2007). The location of the shock depends on 
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FIG. 16. — Contour plots of Nu versus V for the last three time bins of the 
2007 time-resolved spectra. 

the ram pressure of the accretion stream and on the luminos- 
ity, L. For L close to the local Eddington luminosity, with in- 
creasing L the shock moves away from the neutron star, where 
the fi-field is weaker. Thus the CRSF moves to lower ener- 
gies, as is indeed observed in most sources. For smaller L, 
however, the ram pressure of the accreting material becomes 
important. As shown by Staubert et al. (2007) for the case 
of Her X-l, once L is below a critical luminosity the height 
of the column decreases as M (and L) increases, because the 
accreted material "squeezes" the accretion mound. In this 
regime, therefore, the cyclotron line energy increases with lu- 
minosity. 4U 1907+09, however, is more like A0535+262, 
with no observed correlation between the CRSF energy and 
total flux. This could imply that it is somewhere between the 
two regimes described above. 

We see strong variations in the CRSF parameters over the 
pulse phase, as well as in the continuum paramters, particu- 
larly Ecutoff an d r. These clearly indicate that with a change 
in viewing angle we see different aspects of the accretion col- 
umn structure. For details on phase-resolved pulsar spectra 
see Coburn et al. (2002). 

6.2. The Pulse Period Evolution 

As seen in §4, 4U 1907+09 has an interesting and as it turns 
out, rather unique pulse period history. Typical XRBs tend to 
maintain a trend of either spin-up or spin-down with short, in- 
tense episodes of the contrary superimposed (Bildsten et al. 
1997). Examples of such behavior include Her X-l, Cen X-3, 
Vela X-l, and X Per. These shorter episodes are often inter- 
preted to be due to short term variations in the mass transfer 
rate from the donor star which give rise to short-term torque 
fluctuations (Ghosh & Lamb 1979). Such changes are always 
accompanied by large magnitude changes in L. The model 
used to describe this behavior is a magnetic coupling between 
the neutron star and its truncated accretion disk which trans- 
fers angular momentum to the neutron star, the accreted mat- 
ter producing a torque and giving rise to non-zero values of 
P (see, e.g., Ghosh, Pethick & Lamb 1977, and references 
therein). 

There are only three sources known to show a very differ- 
ent behavior: The LMXB GX 1+4, the ultracompact LMXB 
4U 1626-67, and 4U 1907+09, all three of which were his- 
torically characterized by nearly constant P over decade long 
periods with no short-term fluctuations but rare dramatic re- 
versals of P (Chakrabarty et al. 1997a,b; Bildsten et al. 1997; 
Fritz et al. 2006). Monitoring of GX 1+4 by the Fermi 
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Gamma-ray Burst Monitor show it continuing in its new trend 
6 , however 4U 1626-67 has shown another torque reversal 
(Camaro-Arranz et al. 2009), returning to a spin up trend af- 
ter 18 years of spin down. Despite the difference of HMXB 
versus LMXB, 4U 1907+09 and 4U 1626-67 have shown 
very similar pulse period behavior, both with \P/P\~ several 
1000 years, neither showing any change in X-ray luminosity 
during reversals of P, and both displaying similar P/P ratios 
for their spin up and spin down periods. 

With the simple magnetic torquing model outlined above, 
both the long distinct episodes of constant P and the mag- 
nitude and sign of P are difficult to explain together (Fritz 
et al. 2006, and references therein). As was observed in 
GX 1+4 (Chakrabarty et al. 1997a), one would expect dra- 
matic changes in luminosity to accompany any torque rever- 
sal. No such correlation has been seen in either 4U 1626-67 
or4U 1907+09. Furthermore, it is difficult to understand why 
the torquing in the spin up regime has a similar magnitude but 
reversed sign than in the spin down regimes. 

Recently, Perna, Bozzo & Stella (2006) presented a new 
model in which torques can change without having to invoke 
retrograde disks or fluctuations in M. This oblique rotator 
model can account for all observed properties of the spin his- 
tory of 4U 1626-67, including the large values of \P/P\ be- 
fore and after the torque reversal and the virtually unchanged 
luminosity of the source. However its prediction that such 
torque reversals must be very rare is contradicted by the sec- 
ond observed reversal of 4U 1907+09 in only a few years. 
Further monitoring of the pulse period will be needed as the 
evolution of the period of 4U 1907+09 may help to improve 
models of the accretion mechanism in XRBs. 

6.3. The Fe lines and the Absorbing Gas: Tracing the 
Accretion Flow 

We have studied the Fe K emission complex in detail, re- 
solving the Fe Ka line, detecting Fe Kf3 for the first time 
in this object, and detecting a strong Fe K edge at 7.11 keV 
which we model employing an Fe abundance 4 times the solar 
value in the absorbing material. 

It is possible that the same material which absorbs the X- 
ray continuum below 2 keV along the line of sight is respon- 
sible for transmitting the fluorescent emission lines. We can 
assume for simplicity that the absorbing gas is situated in 
an optically-thin, uniform-density, spherical shell of material 
surrounding the central X-ray source. We can thus use Eqn. 
4 of Yaqoob et al. (2001; based on Krolik & Kallman 1987), 
which yields an estimate for the equivalent width, EW ca i c , of 
the Fe Ka line given the photon index of the underlying con- 
tinuum, the column density Nh, assumed Fe abundances, and 
a covering fraction / c of the sky as seen from the central X-ray 
source. 

With the abundances of Wilms et al. (2001), solar abun- 
dance for Fe means 2.69 x 10~ 5 Fe atoms for every H atom. 
We see that for the material obscuring 4U 1907+09, the Fe 
abundance is ^4 times this value. Using the best-fit values 
of T (1.15) and A^ H (2.18 x 10 22 cm" 2 ) from the 2007 time- 
averaged FDCut fit, we find EW ca , c = 47 eV / c . With an ob- 
served EW of 47±4eV, we can conclude that the Fe Ka line 
is consistent with being transmitted by the absorbing material 
if the covering fraction is close to unity. It is thus likely that 
the absorbing material along the line of sight is in the stellar 

6 http://gammaray.msfc.nasa.gov/gbm/science/pulsars/ 



wind in the vicinity of the pulsar. The high covering factor im- 
plied for the Fe line emitter, the lack of continuum emission 
observed below 1 keV in the time-averaged spectrum, and the 
absorbing/Fe-line emitting material being full-covering are all 
consistent with a wind accreting geometry, as opposed to a 
disk accreting geometry. 

The Fe Ka line energy, the K/3/Ka intensity ratio, and the 
properties of the absorber (e.g., edge consistent with neu- 
tral Fe K) all suggest that the bulk of the absorbing/Fe-line 
transmitting material is neutral or at most moderately-ionized. 
Theoretical calculations based on Kallman & McCray (1982), 
however, indicate that the stellar wind should be fully ionized 
at the radius of the neutron star. According to Kallman & 
McCray iron becomes fully ionized above log(Q = 3 where 
£= L/nR 2 . Using the stellar wind parameters from Cox et al. 
(2005) we can calculate the ionization parameter, log((,) = 3.9, 
indicating full ionization. We calculate the distance from the 
companion star at which the ionization of the wind is con- 
sistent with the observed emission features to be roughly 16 
times the orbital radius of the pulsar. The existence of neu- 
tral, line-emitting iron close to the neutron star in the midst of 
a highly ionized wind could be possible if the wind from the 
companion star is clumpy rather than homogeneous. Sako et 
al. (1999) used the clumpy wind model to reconcile the pres- 
ence of flourescent lines in the spectrum of the HMXB Vela 
X-l by considering an inhomogeneous wind with cool, dense 
clumps of neutral or near-neutral material. We will show in 
§6.4 using time variability arguments that this model is likely 
for the accretion in 4U 1907+09. 

6.4. Variability of the Continuum Emission and Absorbing 
Components 

Time-resolved spectral modeling shows that variations in 
the power-law component are observed in the 1-40 keV con- 
tinuum variability. The multi-band lightcurves show corre- 
lated flux variations up to 10 keV, in agreement with this 
notion. Given the low column density, variations due to the 
absorbing material cannot account for continuum variability 
above 3 keV, particularly in the prolonged dip in the 2006 ob- 
servation and the flares in the 2007 observation; though as 
noted earlier, the change in T coinciding with the dip at the 
end of the 2007 observation may be due to incomplete mod- 
elling of a partial-covering absorber as indicated by Fig. 1131 

The X-ray continuum in wind-accreting pulsars is generally 
thought to originate in the accretion column as the accretion 
flow is channeled by the neutron star's magnetic fields onto its 
poles and then is decelerated in a radiative shock and settles 
onto the surface of the neutron star. Recent models attribute 
the X-ray emission to inverse Comptonization of soft seed 
photons (originating in the accretion mound and the shock 
front) as the photons interact with the compressing gas (e.g., 
Becker & Wolff 2007). Variations in F could thus indicate 
changes in the optical depth, temperature, or geometry and 
therefore the average number of upscatterings a typical pho- 
ton undergoes. 

Finally, the observation of a doubling in Nh during the fi- 
nal ~30 ks of the 2007 observation further indicates that the 
absorbing material is indeed clumpy. According to recent 
calculations, winds of O and B stars show strong clumping, 
due to instabilities in the line-driven acceleration mechanism 
(Dessart & Owocki, 2005). Oskinova et al. (2007) confirmed 
this result, showing that clumping is indeed necessary to rec- 
oncile measured line strengths with observed mass loss rates 
of the companion star. Observationally, in the HMXB Vela 
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X-l, a system similar to 4U 1907+09, Sako et al. (1999) 
showed that a clumpy wind is required to explain the lines 
seen in the Chandra gratings spectrum. This result was con- 
firmed in a statistical analysis of the Vela X- 1 flaring behavior 
(Fuerst et al., 2008), who found typical clump masses around 
2.3xl0 20 g. 

Flares seen in 4U 1907+09 on the scale of ^lhr as in 
our 2007 observation are consistent with flaring behavior de- 
scribed by Owocki et al. (2009) and Walter & Zurita-Heras 
(2007). Based on the clumpy wind models, Walter & Zurita- 
Heras calculated clump masses from the brightness and dura- 
tion of typical flares. Their results are on the order of 10 21 g 
with A^h values for the clumps around 10 22 cm" 2 . Using sim- 
ilar methods we obtain clump masses of around 1.5 x 10 20 g 
and a value for A^h also around 10 22 cm" 2 . Additionally, this 
value for A^h is consistent with the dip in the 2007 observation 
where the increase in A^h is about 2 x 10 22 cm~ 2 which could 
be caused by a clump passing through the line of sight. 

From the color-color diagrams we can see that the dip in 
the 2007 observation and the overall dimness of the 2006 ob- 
servation are consistent with a partial covering absorber in the 
line of sight. This could be caused by the clumps in the stel- 
lar wind partially obscuring the object. The prolonged dip in 
the 2006 observation does not fit this explanation, however it 
can still be explained by a clumpy wind model. Since clumpy 
winds are expected to have regions where the wind density 
is significantly below its average density this would explain 
dips due to sudden decreases in the mass accretion rate as in 
2006 and previously observed in 4U 1907+09 by in't Zand et 
al. (1997). The clumpy wind model therefore explains the ob- 
served dipping and flaring behavior as well as the presence of 
near-neutral iron in the otherwise highly ionized stellar wind. 

7. CONCLUSIONS 

The use of the Suzaku observatory for our two observations 
has lead to several advantages over previous observations of 
4U 1907+09 in that the broadband coverage of the XIS+PIN 
has allowed us to probe in depth the spectral variability in this 
source. The XIS offers CCD-resolution and a high effective 
area to study Fe K bandpass features, plus higher sensitivity 
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below 3 keV than previous investigations to accurately quan- 
tify absorption. The PIN has very low background covering 
the 1 3-^-0 keV band, allowing us to analyse characteristics of 
theCRSFat 19keV. 

Simultaneous multi-band dips and flares are likely caused 
by changes in accretion rate possibly due to inhomogeneous, 
clumpy winds. Dips in the soft band could also be explained 
by clumpy winds when a clump passes through the line of 
sight. 

We have studied the Fe K a line with unprecedented preci- 
sion and the Fe K/3 line for the first time. We have determined 
that the iron-line emitting material is very close to the pulsar 
and therefore could also be due to cooler, denser clumps in an 
otherwise fully ionized stellar wind. 

The pulse period has undergone dramatic changes recently 
and we are able to fill in the gap between Fritz et al. (2006) 
and inam et al. (2009) where the second torque reversal took 
place. The torque reversals cannot be fully explaned by any 
current model. 

Finally, we have observed the CRSF at ^19keV in our 
time-averaged, phase-resolved and time-resolved spectral fit- 
ting. Lack of variation of the cyclotron line energy with 
flux places 4U 1907+09 in the same regime as A0535+262, 
somewhere in between the large and small luminosity regimes 
where perhaps the effect of the ram pressure is just enough to 
keep the shock front of the accretion column at a fairly con- 
stant height. 
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